/*=============================================================================
* Copyright (C) 2020 All rights reserved.
*
*   Filename:		propagate.hpp
*
*   Author: Wang Zhecheng - wangzhecheng@yeah.net
*
*   Date: 2020-03-24 04:18
*
*   Last Editors: Wang Zhecheng - wangzhecheng@yeah.net
*
*   Last modified:	2020-08-04 01:26
*
*   Description:
*
=============================================================================*/
#ifndef PROPAGATE_H
#define PROPAGATE_H

#include "../geometry/shape/subdivide.hpp"

namespace geoxfem {
namespace xfem {

template <typename Model>
void initialCrack(const std::string& shape_name,
                  const VectorXR&    coord0,
                  const VectorXR&    coord1,
                  Model&             model) {
  using Crack             = typename Model::Crack;
  using Shape             = typename Model::Shape;
  using point_type        = typename traits::point_type<Model>::type;
  constexpr size_t degree = Shape::Degree::value;

  VectorXR vec = coord1 - coord0;

  std::vector<point_type> plist, _plist;
  plist.emplace_back(point_type{coord0});

  for(size_t i = 0; i != degree; ++i) {
    point_type pi(coord0 + static_cast<_real>((i + 1.0) / degree) * vec);
    plist.emplace_back(pi);
  }
  /* 1. create new crack in model */
  Crack* crack = model.append(Crack{model.getMesh()});
  /* 2. add points and segment into model */
  for(auto& point : plist) {
    point_type p = model.append(point);
    _plist.emplace_back(p);
  }

  std::unique_ptr<Shape> new_segment =
    Shape::getFactory().create(shape_name, _plist);
  Shape* p_segment = model.append(new_segment, crack);
  /* 3. sub divide main mesh */
  algorithm::SubDivide(p_segment->toSegment(), *model.getMesh());
  return;
}

template <typename Tip, typename Model>
void propagate(const _real& delta_theta,
               const _real& length,
               Model&       model,
               Tip&         tip) {
  using DimensionType        = typename Model::DimensionType;
  using degree               = typename Model::Degree;
  using Shape                = geometry::ShapeInterface<DimensionType, degree>;
  using point_type           = typename traits::point_type<Shape>::type;
  constexpr size_t dimension = traits::dimension<DimensionType>::value;

  /* get old tip point and new tip point */
  point_type old_tip = tip.getTip();
  _real      theta   = tip.getTangentDirection() + delta_theta;
  VectorXR   delta(dimension);
  delta << std::cos(theta) * length, std::sin(theta);
  point_type new_tip = model.append(point_type(old_tip + point_type(delta)));
  std::unique_ptr<Shape> sp_new_vertex =
    Shape::getFactory().create("Vertex", {new_tip});
  /* generat segment */
  model.append(sp_new_vertex, tip.getParent());
  std::unique_ptr<Shape> new_segment =
    Shape::getFactory().create("segment", {old_tip, new_tip});
  Shape* p_segment = model.append(new_segment, tip.getParent());
  /* devide main mesh */
  algorithm::SubDivide(p_segment->toSegment(), *model.getMesh());

  return;
}

}  // namespace xfem
}  // namespace geoxfem

#endif /* PROPAGATE_H */
